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We investigate transition rates between different spin configurations for S > 1 spins weakly 
coupled to a d-dimensional phonon bath. This study is motivated by understanding observed mag- 
netization relaxation as a function of temperature in diverse magnetic systems such as arrays of 
magnetic nanoparticles and magnetic molecules. We assume that the magnetization of the spin 
system relaxes through consecutive emission or absorption of a single phonon. From a weak, lin- 
ear spin-phonon coupling Hamiltonian, we derive transition rates that would be used to examine 
dynamic properties of the system in kinetic Monte Carlo simulations. Although the derived phonon- 
assisted transition rates satisfy detailed balance, in the case of two and three dimensional phonon 
baths, transitions between degenerate states are not allowed. Thus, if there are no alternative paths 
along which the spin system can relax, the relaxation time diverges. Otherwise, the system finds 
other paths, which leads to an increase in the relaxation time and energy barrier. However, when 
higher-order phonon processes are included in the transition rates, it is found that the system can 
reach the states which were inaccessible due to the forbidden transitions. As a result, the system 
recovers some of the dynamic properties obtained using the Glauber transition rate. 

PACS numbers: 75.60. Jk,02.70.Tt,61.20.Lc 



I. INTRODUCTION 

In many physical and chemical nanoscale systems rang- 
ing from semiconductor quantum dots to arrays of mag- 
netic nanoparticles or nanoscale magnetic molecules, dy- 
namic properties play a crucial role in understanding the 
underlying physics and in designing systems of interest 
for practical applications. For example, the time evo- 
lution of quantum systems into decoherence needs to be 
fully understood in various local environments in order to 
build scalable quantum computers. It is also important 
to investigate the spin-lattice relaxation time Ti and the 
magnetization relaxation time for recently synthesized 
nanoscale magnetic systems to use them as information 
storage devices. 

To study the dynamic properties of the systems dis- 
cussed above, it is common to consider interactions of 
the systems with their environment. The environment is 
typically described as a heat bath, which has a much 
shorter relaxation time than the systems. Depending 
on the physical quantities to be calculated, one has to 
choose an appropriate bath. To understand decoherence 
mechanisms in quantum dots at low temperatures, one 
often considers an interaction between electron spins and 
a spin bath consisting of a large number of 5 — 1/2 
spins that mimic nuclear spins i^i^ To estimate the spin- 
lattice relaxation time in quantum dots, one should take 
into account electronic spins coupled to a phonon bath.'^ 
In nanoscale magnetic systems such as arrays of mag- 
netic nanoparticles^i^ single-molecule magnets^iii^ and 
single-chain magnets )2ii2iiii the effect of nuclear spins is 
minimal, so a phonon bath becomes more relevant than 
a spin bath. For a single-molecule magnet embedded 
in a three-dimensional lattice, transition rates between 
different spin configurations have been derived from cou- 
pling to the lattice via magneto-elastic coupling i^iiiSiii^ 



The magnetization relaxation time was estimated using 
these transition rates and quantum tunneling rates, be- 
ing in good agreement with experimental dataii For a 
nearest-neighbor interacting ferromagnetic Ising system, 
transition rates were derived from a weak, linear cou- 
pling of the system to a one-, two-, or three-dimensional 
phonon bath i^'^i^'^i^^ Using these phonon-assisted tran- 
sition rates, kinetic Monte Carlo simulations were per- 
formed to measure the lifetime of the metastable state or 
magnetization relaxation time at low temperatures. The 
Monte Carlo simulations revealed that the dynamic prop- 
erties obtained using the phonon-assisted transition rates 
greatly differ from those using other transition rates, 
such as Glauber— or Metropolisj^i It is known that the 
Glauber transition rate can be derived from a coupling 
of a spin system to a fermionic bathJ^ Recently, it was 
shown that the soft Glauber transition rat©i2, requires 
different interpretation in the form of the lifetime of the 
metastable state at low temperatures, although the en- 
ergy barrier to reach equilibrium is the same as that for 
the standard hard Glauber transition rate.— Therefore, 
selection of a proper bath and relevant transition rates 
is critical in understanding dynamic properties. The 
transition rates derived by Park et ali^i^^^ (two- and 
three-dimensional baths) and other groupsii^ii^ (three- 
dimensional bath) using coupling to a phonon bath, share 
a common feature that the rates become zero for degen- 
erate states. Very recently, the derived phonon-assisted 
transition rates were used to examine the nanostructure 
of field-driven solid-on-solid interfacesi^ It was found 
that the phonon-assisted rates provide significant differ- 
ences from other types of transition rates, such as the 
Glauber dynamics. 

In this study, targeting arrays of weakly interacting 
magnetic nanoparticles, single-molecule magnets, and 
single-chain magnets, we generalize the formalism used 
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for an Ising system in Refs fTMT^flsl to a 5 > 1 model on 
a lattice. Each spin in the model interacts via its near- 
est neighbors and has easy-axis single-ion anisotropy. We 
assume that all spins in the model are weakly coupled to 
a phonon bath in d — 1,2, or 3 dimensions. Consider- 
ing that spin relaxation occurs through first-order one- 
phonon emission or absorption processes, we derive tran- 
sition rates for one-, two-, and three-dimensional phonon 
baths. In the cases of d = 2 and 3, some transitions or 
relaxation paths are inaccessible because transition rates 
between degenerate states vanish. This results in increas- 
ing magnetization relaxation time and energy barrier to 
be overcome. However, when higher-order phonon pro- 
cesses are included, other relaxation paths are opened up, 
leading to shortening of the relaxation time compared to 
the first-order phonon processes. The formalism for the 
phonon- assisted transition rates is presented in Sec. II. 
The consequences of using the derived transition rates 
are discussed in the context of kinetic Monte Carlo simu- 
lations in Sec. III. Higher-order processes and their effects 
on dynamic properties are presented in Sec. IV. The con- 
clusion follows in Sec.V. 
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FIG. 1; Schematic diagram of (a) one-phonon absorption pro- 
cess for rotating the second spin from m2 to m2 — 1, and (b) 
one-phonon emission process for rotating the third spin from 
ma to ms — 1 for a five-spin S — 1 system. The wavy lines 
denote the phonons emitted or absorbed, (a) The energy 
difference {Ek — Em) equals 2J + D + H that is positive, 
and (b) {Ek — Em) equals J + D + H that is negative when 
-{2.7 + D) <H < -{.J + D). 



II. FORMALISM FOR PHONON-ASSISTED 
TRANSITION RATES 

Although the current formalism can be applied to more 
general cases, we start with the following Hamiltonian for 
Ns spins {S > 1) on a lattice. 



{r,3) 1=1 ^=l 



(1) 



where Sf is the z component of spin operator S at site i, 
J(> 0) is an exchange coupling constant between nearest- 
neighboring spins, and the first summation runs over all 
nearest neighbor pairs. A positive value of J implies 
ferromagnetically coupled spins. D{> 0) is a single-ion 
magnetic anisotropy parameter determined by the spin- 
orbit coupling. A positive value of D indicates that the 
magnetic easy axis of an individual spin is along the ±z 
axis. Notice that our convention on D differs from other 
works. H is an external magnetic field applied to the spin 
system. We call the spin Hamiltonian Eq. ([1]) a S* > I 
model. The eigenstates \ni) of the spin Hamiltonian arc 
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where \mi) is the eigenstate of and — —S,—S + 
1, ...,S - 1,S. For 5 = 1 and D < 0, the spin Hamil- 
tonian Eq. ^ is known as the Blume-Capel model^^i^ 
which was introduced to describe features of the phase 
diagram of He^-He^ mixtures as well as to understand 
a phase transition in UO2. The S > I model can be 
applied to the following magnetic systems: arrays of 
weakly-interacting magnetic nanoparticles^^i^ nanoscale 
single- molecule magnets such as Mni2 and Fegf^i^ a 



Mn(III)2Ni(II) single-chain magnetjiSiiii and a Co fer- 
rimagnetic compoundi^ 

To understand spin relaxation in the S > 1 model, we 
assume that the spin system is weakly linearly coupled 
to a surrounding phonon bath in d dimensions. Since the 
phonon bath has much shorter relaxation time than the 
spin system, it is assumed that each spin is independently 
coupled to the bath. Then spin relaxation occurs through 
consecutive emission or absorption of a single phonon 
with energy that equals the cost of rotating a single spin 
from rrii to rrii ± 1 in a single transition while keeping the 
rest of the spins fixed. Henceforth we call these rotations 
first-order one-phonon processes. 

A phonon bath is described as a collection of simple 
harmonic oscillators so the phonon Hamiltonian is writ- 
ten as 
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where q is the phonon wave vector, uj^ is the angular fre- 
quency of a harmonic oscillator, ct and are creation 
and annihilation operators of a phonon with wave vector 
q, and c is the sound velocity in the lattice. The following 
spin-phonon coupling Hamiltonian Tigp-ph has the sim- 
plest form that takes into account all possible first-order 
one-phonon processes. 



N 



sp — ph 
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where A is a coupling constant, N is the number of unit 
cells associated with the phonon bath, M is the mass 
of the particle in the unit cell, and Sf' are the raising 
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and lowering spin operators for site j. Here the polariza- 
tion of the phonons is not considered for simplicity. The 
magneto-elastic coupling theory^^ suggests that the spin- 
phonon coupling must be proportional to a linear strain 
tensor, ea'a = Vq'Uq, where Ua is the a component of the 
displacement vector u and a', a € {x, y, z}. The displace- 
ment vector can be expressed in terms of ct and c^, and 
a Fourier transform is carried out on ea'a- This explains 
the dependence of the prefactor of 7isp-ph on the wave 
vector q. Equation ((4]) contains the minimum number of 
terms required to rotate the spin vectors via onc-phonon 
emission or absorption processes. Notice that the spin 
system includes a nearest-neighbor exchange interaction, 
in contrast to works reported by other groups^S. Due to 



the exchange interaction, a transition from mi to mi — 1 
does not uniquely determine the sign of the energy differ- 
ence Ek — Em between the two states. For example, when 
mi = 1, Ek — Em is positive for Fig. [IJa), while Ek — Em 
is negative for Fig. [TJb). The same rule is applied to 
transition from m; to rn; -f 1. 

Using Fermi's golden rule, within perturbation the- 
ory, we calculate the transition rate from state |m) = 
mi) (g) \m2) (8) • • • (8) ■ • • I'itT'n) to state |fc) = 
mi) \m2) • • • (8) I'lTT'i) ■ • • If^N), where m[ = 
TO; ± 1. States |m) and |fc) differ by the rotation of a 
single spin (m; m; ± 1) at site I. We first consider the 
transition rate Wkm from |m) to |fc) for emission of one 
phonon with energy Hco^, as illustrated in Fig. [T]^b) . 



Wkm = l,k\nsp-ph\nq,m)\'^ pphS{huJi^- {Em - Ek)), 
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where pph is the phonon density of states, is the 
eigenvalue of the phonon number operator, and ctjng-) = 

^JTl^'+\\n^+l) is used. Here Ek and Em are the energies 
of states I A:) and |m) calculated from the spin Hamilto- 
nian, Eq. ([1]), and the energy difference Ai? is given by 



AE{mi ^ m; ± 1) 



= Ek 



TD{2mi±l)T(m 



where the sum runs over nearest neighbors of site I. Using 
the Bose-Einstein distribution function, one knows that 
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where fc^ is the Boltzmann constant, T is the temper- 
ature, and (3 — l/{kBT). Assuming that the bath re- 
laxes much faster than the spin system, we integrate over 
all degrees of freedom of the bath and convert J2q i^i^o 
[(A^a'')/(27r)'^] / d'^q, where a is the lattice spacing. Then 
the transition rate from state |m) to |fc) becomes 
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{Em — Ek)'^ 



N 



1 



{S + mi){S -mi + l)5ki^mi-i 
+ {S - mi){S + mi + l)4,_„i,+i, 
2tt {d = 3), 2 (d = 1,2), 



Em - Ek M 

(10) 
(11) 



where is a mass density associated with the bath and 
ki is the quantum number of the ^th spin for state |fc). 



The transition rate for absorption analogously becomes 
A2 TV {Ek - Em)'' 



Wkm = 



Ek - Em >(I)2) 



Henceforth, we refer to the derived transition rates, 
Eqs. ([9]) and (fT2|) . as phonon-assisted transition rates. 



III. CONSEQUENCES OF PHONON-ASSISTED 
TRANSITION RATES 

The derived phonon-assisted transition rates are gen- 
eralized forms of those for the Ising model discussed in 
RefsEHlll According to Eqs. Q and the transi- 
tion rates are highest for mj = and lowest for m; — ±5* 
for large S. A transition rate^ similar to the derived 
rates was obtained for an isolated spin cluster embedded 
in a lattice instead of interacting spin clusters. In this 
formalism, the single-ion anisotropy parameter D corre- 
sponds to the coupling constant A in Eq. Assuming 
that X ^ D, we estimate the magnitude of the prefactor 
of the c? = 3 phonon-assisted transition rates, Eqs. (O 
and (fT2|) . for example, for the nanoscale single- molecule 
magnet Mni2. Using measured parameter values such 
as S' = 10, mi = 10, J = 0.01 K, 7] = 1.83 x 10^ kg/m^ 
[Ref l26j ■ and c = 1.45 x 10'^ m/s^Z we find the prefactor to 
be 0.00041 s"^ Here we use J = 0.01 K due to the large 
intermolecular separation in this system, although it was 
not directly measured. To associate the derived transi- 
tion rates with magnetization relaxation times for various 
nanoscale systems, one needs to solve a master equation 
including the derived transition rates or perform kinetic 
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FIG. 2: Transition rates W vs AE/T for the Glauber and 
phonon-assisted transition rates in the case of d = 1,2, 
and 3 dimensional phonon baths, computed at temperature 
T — l.OJ/fes. The prefactor of the phonon-assisted transition 
rates was not included. 



Monte Carlo simulations with the rates. Hereafter we 
focus on the latter approach. 

In Fig. [2] the phonon-assisted transition rates for d = 
1,2, and 3 dimensional phonon baths arc shown as func- 
tions of AE/T at temperature T — l.OJ and com- 
pared to the Glauber transition rate, l/(e^(^'=~^'"' + 1). 
The main difference between the phonon-assisted and 
Glauber transition rates originates from the nature of 
the bath coupled to the spin system. When degenerate 
states are involved in transitions, this difference becomes 
prominent. In the d ~ 2 and 3 phonon-assisted transition 
rates, transitions between degenerate states are forbidden 
because Wkm = when AE = 0. In the case of c? = 1, 
the rate does not vanish for AE = 0, but rather de- 
creases monotonously with AE/T like the Glauber tran- 
sition rate. Ramifications of the forbidden transitions 
on dynamic and equilibrium properties are discussed in 
Monte Carlo simulations. For simplicity, we consider a 
ferromagnetic S — 1 model [the Blume-Capel model with 
_D > in the spin Hamiltonian, Eq. ([T])] on a L x L square 
lattice with the d = 2 phonon-assisted transition rates 
unless specified otherwise. 



A. Dynamic properties 





H/J 


\ +1 

o\ 







D/J 


/ -1 





FIG. 3: Phase diagram of the Blume-Capel model with D > 
in the spin Hamiltonian ([T| at zero temperature.— The equi- 
librium spin configurations are shown in the three different 
regions. 
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FIG. 4: Schematic diagram of relaxation of magnetization for 
a ferromagnetic S = 1 Blume-Capel model at H — — 3.25J 
and D = 0.75 J, via first-order one-phonon processes, mi 
m; ± 1. State (A) denotes all spins aligned along the +z axis, 
and state (B) a single m; = spin in the sea of rrii = -1-1 spins. 
The boxes represent rotated spins from the +z axis in the ini- 
tial state. State (C) represents two nearest-neighbor m; = 
spins in the sea of rrii = -1-1 spins, while state (G) represents 
two mi — spins that are not nearest neighbors. The numbers 
next to the arrows indicate the energy difference AE defined 
by Eq. With the d = 2 phonon-assisted transition rates, 
the spin system relaxes through (A)^(B)— >(C)^(E) rather 
than (A)— >(B)^(D)^(E). The thick arrows represent the 
most probable path for relaxation with the phonon-assisted 
rates. The thin arrows denote the most probable path with 
the Glauber transition rate. 



We investigate the effects of the forbidden transi- 
tions on the lifetime of the metastable state for the 
5 = 1 Blume-Capel model below the critical tempera- 
ture. The equilibrium spin configurations for the model 
are shown in different regions in FiglpI The nucle- 
ation and mctastability for the model were studied for 
-4 J < D < -3 J and -J < H < -{U + D) in Ref,27. 
The critical temperature Tc increases as D/J increases 
at H=0. The value of Tc at D = was calculated us- 
ing different methods. Bethe-lattice approximation gave 
rise to 2.065 J/ and the effective field theory sug- 
gested 1.952 J/fcs;^ Expanded Bethe-Peierls approxima- 



tion produced 1.915J//cb>'^*' and Monte Carlo simula- 
tions suggested that the critical temperature at H=0 is 
1.6950J//CB at L> = and 2.1855J/A:i3 at D = 5J^ 
Suppose that all spins are initially aligned along the +z 
axis. When an external magnetic field is applied along 
the —z axis, the initial state becomes metastable. When 
— (4J + D) < H < 0, at low temperatures, the spin sys- 
tem relaxes toward the stable state (all spins along the 
—z axis) via creating a single critical droplet consisting 
of connected = — 1 spins. (In the regime studied in 
Ref,'27' multiple critical droplets are formed.) To give a 
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FIG. 6: A checkerboard state for a S' = 1 SDin svstem on a 



FIG. 5: Transitions forbidden by first-order one-phonon pro- 
cesses for a S' = 1 spin system on a square L x L lattice with 
D = 4J and // = in which Wkm ~ 0. (a) Initial state of 
— 1 with nearest neighbors of -f 1. (b) Initial state of 4-1 with 
nearest neighbors of —1. 



specific example, we consider H — — 3.25J, D — 0.75J, 
and T < 0.02J/fcs- Figure 2] illustrates a few possible 
relaxation paths from the metastable state. The ini- 
tial state (A) can first relax to state (B), which repre- 
sents a single mi — spin in the sea of mi = +1 spins. 
Then state (B) can relax to one of the states (C), (D), 
or (G), where state (C) denotes a single mi = —1 spin 
in the sea of to; = +1 spins and state (D) [state (G)] 
two nearest-neighbor mi = spins [two m; = spins 
that are not nearest neighbors] in the sea of m; = +1 
spins. Transitions between states (B) and (D) are for- 
bidden because AE = ±(4J - D + H) = 0. Thus, the 
system can relax from state (B) to either state (C) or 
(G). The most likely path among the alternative ones 
is (B)^(C)^(E), as indicated by the thick arrows in 
Fig. [31 Thus, for the d = 2 phonon-assisted transition 
rate, the critical droplet is state (C). As a result, the en- 
ergy barrier to be overcome in order to reach the stable 
state, is 2.0J, while it is 1.5J for the Glauber transition 
rate. Accordingly, the relaxation time becomes longer 
than that for the Glauber transition rate. In an Ising 
system coupled to a phonon bath, a similar behavior has 
been found at H — —2JJ^^^^^ If there are no alternative 
paths, the relaxation time diverges. For a two-electron 
state in a quantum dot, the measured spin-lattice relax- 
ation time was observed to diverge at a particular mag- 
netic field where the triplet and the singlet states became 
degenerate 1^ 



B. Equilibrium properties 

The d — 1,2, and 3 phonon-assisted transition rates 
satisfy detailed balance. Thus, equilibrium properties ob- 
tained using the derived transition rates must agree with 
those obtained using different transition rates in Monte 
Carlo simulations. However, there is a caveat in this 
statement because of forbidden transitions between de- 
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FIG. 7: Cumulant Ul [Eq. (|13p ] vs temperature calculated 
from Monte Carlo simulations using (top) the Glauber and 
(bottom) the d — 2 phonon-assisted transition rates for S — 1 
spins on a LxL square lattice at D = 4 J and H = 0. Random 
initial states and periodic boundary conditions were used. For 
the fixed point of the cumulant and critical temperature, see 
the text. 



generate states for the d = 2 and 3 phonon-assisted rates. 
In some cases the forbidden transitions would prevent the 
spin system from relaxing to the equilibrium state if we 
start with a particular initial state. As an example, we 
consider D — AJ and H ~ 0. As illustrated in Fig. O 
transitions between mi — +1 {mi — —1) and m; = 
with the sum of the nearest neighbors fixed as —4 (-1-4) 
are not allowed because AE = —4,/ ± H + D ^ 0. So 
if we started with a checkerboard initial state as shown 
in Fig. [6] in Monte Carlo simulations, the system would 
stay indefinitely at the initial state because Wkm — 
for any possible single-spin rotations. However, we have 
confirmed that the system reaches equilibrium if we start 
with a random initial state or a state slightly modified 
from the checkerboard pattern. (Even one defect site in 
the perfect checkerboard state is sufficient.) A similar 
feature was reported in the time evolution of field-driven 
solid-on-solid interfaces using the d — 2 phonon-assisted 
transition ratci^ 

Monte Carlo simulations were performed for D = 4J 
and H = with periodic boundary conditions and a ran- 
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dom initial state. L = 20,30,40,50, and 60 were con- 
sidered at several different temperatures. The average 
absolute magnetization per site and cumulant Ul 

were calculated. 
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where m is the magnetization per site and M is the 
total number of Monte Carlo steps. It was found that 
(|m|) (not shown) and Ul (Fig. [7]) as functions of tem- 
perature agree with those obtained using the Glauber 
transition rate for the different system sizes. For the 
Glauber transition rate, the values of Ul for the differ- 
ent system sizes intersect with one another at [/£=0.611± 
0.008 and rc=2.052±0.005 J/ks, while for the phonon- 
assisted transition rate we have [7^=0. 611± 0.009 and 
rc=2.052±0.006 J/kB- So both transition rates give 
rise to the same critical temperature. This critical tem- 
perature is slightly lower than that reported in Reff3ll. 
2.13J/A:b. 
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FIG. 8: Schematic diagram of relaxation of magnetization for 
a. S = 1 spin system a.t H = — 3.25J and D — 0.75J when 
both first-order (mi — > m; ± 1) and second-order one-phonon 
processes {mi mi it 2) are considered. The thick solid and 
dashed arrows represent two highly probable paths for relax- 
ation. Compared to Fig. |4l the spin system can relax from 
state (A) to state (E) through (B), or to state (F) through 
(D). The numbers next to the arrows indicate the energy dif- 
ference AE defined by Eqs. JT]) and (fTf)! . 



IV. HIGHER-ORDER PROCESSES 

We have, so far, discussed spin relaxation caused by 
first-order one-phonon processes, mi — > m; ± 1. However, 
higher-order processes such as multi-phonon processes 
and second-order one-phonon processes (m; — > m; ± 2), 
can also contribute to the spin relaxation. Their contri- 
butions become significant especially when encountering 
transitions forbidden by the first-order one-phonon pro- 
cesses. In this case, the spin system would find less costly 
paths via the higher-order processes than paths directed 
by the first-order processes. Thus, the relaxation time 
becomes shortened and the system may recover the same 
relaxation time or energy barrier as the Glauber tran- 
sition rate. The system starting with the checkerboard 
state (Fig. [6]) can be also relaxed to equilibrium. 

Recently, multi-phonon processes such as Raman pro- 
cesses and two-phonon processes were considered in a 
spin-phonon relaxation rate for rigid atomic clusters, and 
it was shown that there are no closed analytical forms 
for the rate>^ In the current study, we focus on second- 
order one-phonon processes as higher-order processes. 
These second-order processes were included as a part of 
the relaxation mechanism for the single-molecule magnet 
Mni2.— The simplest form of the spin-phonon coupling 
Hamiltonian for these processes is given by 



N 

(2nd) _ \'\^\^ 

sp-ph - ^ 2^2^ 



2NMlu. 



+ (57fct.+ (5+)%~-F(57)%H, (14) 

where A' is a coupling constant and |A'| <^ |A|. Following 
the method used in Sec. II, we obtain the transition rate 



^fem '^^ from state \m) to |fc), where these two states differ 



by a single spin rotation at site = to; ± 2. 



(2nd) _ {XT N' 



(15) 



N' ^{S + mi){S -mi + l){S + mi- 1){S - mi + 2)4„m,-2 
+ (5 - mi){S + mi + l){S-mi- \){S + m/ + 2)4,,m,+2(16) 

where 7 is defined in Eq. (Ilip . The energy difference 
Ek - Em = A£; is 

AE{mi ^ mi±2) = =f2J ^ mf ^'iD{mi±l) 

k^nn{l) 

T2H, (17) 

where the sum runs over nearest neighbors of site I. This 
formula is applied to both emission and absorption pro- 
cesses. 

When the second-order transition rate W^^^^'' is in- 
cluded in the calculation of the lifetime of the metastable 
state, relaxation scenarios are greatly modified as illus- 
trated in Fig. [51 With the same parameter values used in 
Sec.III.A (D = 0.75J and H = -3.25J), the system can 
now relax through transitions as indicated by the thick 
dashed arrows [(A)^(D)^(F)] or by the thick solid ar- 
rows [(A)-^(B)^(E)] in Fig. [51 None of these transi- 
tions involve degenerate states. The critical droplet for 
the first relaxation route is state (D), while that for the 
second route is state (B). In both relaxation paths, the 
energy barrier is 1.5 J, which is the same as that for the 
Glauber transition rate. In the case of equilibrium Monte 
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Carlo simulations, the second-order processes allow the 
system to relax via alternative second-order transitions 
with lower energy cost but AE ^ 0. So the checkerboard 
initial state (Fig. [S]) can reach equilibrium for D — \J 
and H — 0. It is confirmed that for L = 20 the equi- 
librium properties computed with addition of wj^'^'^'^^ to 
agree with those obtained using the Glauber tran- 
sition rate. 



V. CONCLUSION 

We have considered the 5* > 1 model weakly cou- 
pled to a one-, two-, or three-dimensional phonon bath 
and derived corresponding transition rates from the spin- 
phonon coupling Hamiltonian. The derived phonon- 
assisted transition rates for two- and three- dimensional 



baths differ from other transition rates in that the for- 
mer rates become zero for degenerate states. This caused 
some transitions to be forbidden by the first-order one- 
phonon assisted transition rates, increasing the magneti- 
zation relaxation time. Using a combination of the first- 
order one-phonon processes with the second-order pro- 
cesses, however, the system found more energy-efficient 
paths to equilibrium, and the relaxation time shortened. 
These results represent a major step toward developing 
physically realistic kinetic Monte Carlo simulations for 
magnetic spin systems. 



Acknowledgments 

The author is grateful to M. A. Novotny, L. Solomon, 
and P. A. Rikvold for discussions. 



* Electronic address: lkyungwha@vt.edu| 

^ W. Zhang, V. V. Dobrovitsky, K. A. Al-Hassaniah, E. Do- 

gotto and B. N. Harmon Phys. Rev. B 74, 205313 (2006). 
^ N. V. Prokof'ev and P.C.E. Stamp, Rep. Prog. Phys. 63, 

669 (2000). 

^ T. Meunier, I. T. Vink, L. H. Willems van Beveren, K.-J. 

Tielrooij, R. Hanson, F. H. L. Koopens, H. P. Tranitz, W. 

Wegscheider, L. P. Kouwenhoven, and L. M. K. Vander- 

sypen, Phys. Rev. Lett. 98, 126601 (2007). 
" P. Zhang, F. Zuo, F.K. Urban, A. Khabari, P. Griffiths, A. 

Hosseini-Tehrani, J. Mag. Mag. Mat. 225, 337 (2001). 
^ Y.I. Kim, D. Kim, C.S. Lee, Physica B-Cond. Matt. 337, 

42 (2003). 

® P. Politi, A. Rettori, F. Hartmann-Boutron, and J. Villain, 
Phys. Rev. Lett. 75, 537 (1995) 

M. N. Leuenberger and D. Loss, Phys. Rev. 61, 1286 
(2000). 

® E. M. Chudnovsky, D. A. Garanin, and R. Schilling, Phys. 

Rev. B 72, 094426 (2005). 
^ A. Caneschi, D. Gatteschi, N. Lalioti, R. Sessoli, L. Sorace, 

V. Tangoulis, and A. Vindigni, Chem-Eur. J. 8, 286 (2002). 

C. Coulon, R. Clerac, L. Lecren, W. Wernsdorfer, and H. 

Miyasaka, Phys. Rev. B 69, 132408 (2004). 
" J. Kishine, T. Watanabe, H. Deguchi, M. Mito, T. Sakai, 

T. Tajiri, M. Yamashita, and H. Miyasaka, Phys. Rev. B 

74, 224419 (2006). 

L. Solomon, M.A. Novotny, S. Caliskan, K. Park, and P.A. 
Rikvold, submitted to Phys. Rev. B. 

M. A. Novotny and K. Park, Computer Phys. Commun. 
147, 737 (2002). 

K. Park and M. A. Novotny, in Computer Simulations 
Studies in Condensed Matter Physics XIV, edited by D.P. 



Landau, D. Lewis, and H.-B. Schiittler, (Springer- Verlag, 
Berlin, 2002), p. 182. 

K. Park, M. A. Novotny, and P. A. Rikvold, Phys. Rev. E 
66, 056101 (2002). 

R. J. Glauber, J. Math. Phys. 4, 294 (1963). 

N. Metropolis, A. Rosenbluth, M. Rosenbluth, A. Teller, 

and E. Teller, J. Chem. Phys. 21, 1087 (1953). 

P.-A. Martin, J. Stat. Phys. 16, 149 (1977). 

P. A. Rikvold and M. Kolesik, J. Phys. A: Math. Gen. 35, 

L117 (2002). 

^° K. Park, P. A. Rikvold, G. M. Buendi'a, and M. A. Novotny, 
Phys. Rev. Lett. 92, 015701 (2004). 

K. Saito, S. Takesue, and S. Miyashita, Phys. Rev. E 61, 
2397 (2000). 

G. M. Buendi'a, P. A. Rikvold, M. Kolesik, K. Park, and 

M. A. Novotny, Phys. Rev. B 76, 045422 (2007). 
23 H.W. Capel, Physica 32, 966 (1966). 

M. Blume, Phys. Rev. 141, 517 (1966). 
2^ E. Callen amd H. B. Callen, Phys. Rev. 139, A455 (1965). 

T. Lis, Acta Crystallogr. B 36, 2042 (1980). 

E.N.M. Cirillo and E. Olivieri, J. Stat. Phys. 83, 473 

(1996). 

2* Y. Tanaka and N. Uryu, J. Phys. Soc. Japan 50, 1140 
(1981). 

23 H. Polat, U. Akinci, and I. Sokmen, Phys. Stat. Sol. B 240, 
189 (2003). 

3° A. Du, Y.Q. Yii, and H.J. Liu, Physica A 320, 387 (2003). 
3^ R. da Silva, N. A. Alves, and J.R. Drugowich de Felfcio, 

Phys. Rev. E 66, 026130 (2002). 
32 C. Calero, E. M. Chudnovsky, and D. A. Garanin, Phys. 

Rev. B 74, 094428 (2006). 



